perm filename X[S,AIL] blob sn#072709 filedate 1973-11-18 generic text, type T, neo UTF8
COMMENT ⊗   VALID 00003 PAGES 
RECORD PAGE   DESCRIPTION
 00001 00001
 00002 00002	DSCR POW, FPOW, LOGS, FLOGS.  BOTH RETURN REALS.
 00007 00003	DSCR EXP,ALOG -- FOR USE BY EXPONENTIATION ROUTINES & WORLD
 00015 ENDMK
⊗;
DSCR POW, FPOW, LOGS, FLOGS.  BOTH RETURN REALS.
SID  CLOBBERS LPSA,TEMP,USER
CAL SAIL
DES CALLS GENERATED BY COMPILER FOR ↑ OPERATOR

FPOW:	REAL←FPOW(INTEGER!EXPONENT,REAL!BASE)
POW:	REAL← POW(INTEGER!EXPONENT,INTEGER!BASE)

LOGS:	REAL← LOGS(REAL!EXPONENT,INTEGER!BASE)
FLOGS:	REAL←FLOGS(REAL!EXPONENT,REAL!BASE)

SPECIAL CASES:
	A↑0 = 1
	0↑B = 0 IF B GEQ 0.
	0↑B = INF. IF B<0  ; MESSAGE PRINTED
	A↑B = (-1)↑B*|A|↑B  IF A<0, B INTEGRAL
	A↑B = REALPART(A↑B) IF A<0, B NOT INTEGRAL ; MESSAGE

MESSAGE IS PRINTED IF OVERFLOW OR UNDERFLOW HAPPENS.
IN THIS CASE, FIXUP IS MADE SO THAT ANSWER IS EITHER 0, +INF, OR
-INF.

⊗
IFN ALWAYS,<	BEGIN	UTILS>


HERE(FPOW)
	SKIPA 	USER,-1(P)	;BASE
HERE(POW)
	FLOAT	USER,-1(P)
FPX:	MOVM	LPSA,-2(P)	;GET ABS(EXPONENT)
	JUMPE	LPSA,EXZERO	;0 EXPONENT
	MOVSI	A,(1.0)		;SET FOR FLOATING	
	JRST	2,@[FEXS]	;CLEAR AR FLAGS

FEXL:	ASH	LPSA,-1		;PREPARE TO LOOK AT NEXT BIT.
	FMPR	USER,USER	;SQUARE BASE
	 JFOV	 FPOWOV		;OVERFLOW/UNDERFLOW
FEXS:	TRZE 	LPSA,1		;COLLECT PRODUCT?
	FMPR	A,USER		;YES
	 JFOV	 FPOWOV		;OVERFLOW?
	JUMPN	LPSA,FEXL	;LOOP UNTIL EXPONENT ZERO.
	SKIPGE	-2(P)		;POSITIVE EXPONENT?
	   JRST	FEXDU1
POWRET: SUB	P,X33
	JRST 	@3(P)

FEXDU1:	MOVM	LPSA,A		;CHECK FOR OVERFLOW POSS.
	CAMGE	LPSA,[XWD 2400,1] ;SMALL NUMBER 
	 JRST	 FPDOV		;CALL UNDERFLOW
	MOVSI	LPSA,(1.0)	;TAKE RECIPROCAL OF ANS.
	FDVRM	LPSA,A
	JRST	POWRET		;AND RETURN IT.

EXZERO:	SKIPN	USER		;0↑0
ZRET:	 TDZA	 A,A		;RETURN 0
	MOVSI	A,(1.0)		;RETURN FLOATING 1
	JRST 	POWRET


FPOWOV:	SKIPN	TEMP,OVPCWD	;IF TRAPS ENABLED, USE EM
	 JSP	 TEMP,.+1	;ELSE GET FLAGS THIS WAY
	TLNE	TEMP,100	;SKIP IF NOT UNDERFLOW
FPDOV:	 MOVNS	 -2(P)		;UNDERFLOW -- CHANGE EXPONENT SIGN.
	MOVE	A,[XWD 400000,1] ;LARGE NEGATIVE NUMBER
	SKIPG	TEMP,-2(P)	;CHECK SIGN OF EXPONENT.
	 MOVEI	 A,0		;NEGATIVE ==> RESULT 0.
	SKIPGE	-1(P)		;CHECK SIGN OF BASE.
	 TRNN	 TEMP,1		;XOR SIGN OF EXPONENT.
	 MOVNS	 A		;MAKE +- LARGE NUMBER
	ERR	<Exponentiation under or overflow>,1
	JRST	POWRET		;RETURN.

HERE(FLOGS)
.FLOGS:	SKIPA	USER,-1(P)	;FLOATING BASE
HERE(LOGS)
.LOGS:	FLOAT	USER,-1(P)	;FLOAT THE BASE
	SKIPN	-2(P)		;IF ZERO EXPONENT,
	 JRST	 EXZERO		;GO TO COMMON CODE.
	MOVM	TEMP,-2(P)	;CHECK TO SEE IF 'FIX' WILL 
	CAMLE	TEMP,C1		;OVERFLOW
	 JRST	 USLGEP		;YES -- GO TO LOG-EXP
	FIX	TEMP,-2(P)	;CHECK TO SEE IF EXPONENT
	FLOAT	LPSA,TEMP	;HAPPENS TO BE AN INTEGER
	CAMN	LPSA,-2(P)	;IF SO, USE LOOPS TO
	 JRST	 [MOVEM TEMP,-2(P) ;BE SURE OF CORRECT SIGN
		  JRST FPX]
USLGEP:	JUMPE	USER,[SKIPGE -2(P) ;IF BASE ZERO, AND EXPT NEG.
			JRST FPDOV ;RETURN LARGE NUMBER
			JRST ZRET] ;ELSE RETURN ZERO.
	PUSH	P,USER		;ARGUMENT TO 'ALOG'
	PUSHJ	P,.LOG		;CALL IT.
	FMPR	A,-2(P)		;MULTIPLY BY EXPONENT
	PUSH	P,A		;ARGUMENT TO 'EXP'
	PUSHJ	P,.EXP		;CALCULATE
	JRST	POWRET		;AND RETURN.

C1:	243777777777		;2↑35 - EPSILON


DSCR EXP,ALOG -- FOR USE BY EXPONENTIATION ROUTINES & WORLD
SID	CLOBBERS LPSA,TEMP,USER
CAL	SAIL
⊗


;FLOATING POINT SINGLE PRECISION EXPONENTIAL FUNCTION
;THE ARGUMENT IS RESTRICTED TO THE FOLLOWING RANGE
;	-88.028<X<88.028
;IF X<-88.028, THE PROGRAM RETURNS ZERO AS THE ANSWER
;IF X<88.028, THE PROGRAM RETURNS X AS THE ANSWER
;THE RANGE OF THE ARGUMENT IS REDUCED AS FOLLOWS:
;EXP(X) = 2**(X*LOG(E)BASE2) = 2**(M+F)
;WHERE M IS AN INTEGER AND F IS A FRACTION
;2**M IS CALCULATED BY ALGEBRAICALLY ADDING M TO THE EXPONENT
;OF THE RESULT OF 2**F. 2**F IS CALCULATED AS

;2**F = 2(0.5+F(A+B*F↑2 - F-C(F↑2 + D)**-1)**-1

;THE ROUTINE HAS THE FOLLOWING CALLING SEQUENCE:
;	PUSH P,ARG
;	PUSHJ	P,EXP
;THE ANSWER IS RETURNED IN ACCUMULATOR A

HERE(EXP$)
.EXP:	PUSH	P,[0]		;ONE WORKING CELL
	PUSH	P,B		;AND ONE SAVED AC
	MOVE	LPSA,-3(P)	;GET ARGUMENT
	MOVM	A,LPSA		;GET ABSF(X)
	CAMG	A, E7		;IS ARGUMENT IN PROPER RANGE?
	JRST	EXP1		;YES, GO TO ALGORITHM
	ERR <EXP: under or overflow>,1
	HRLOI	A, 377777	;GET LARGEST FLOATING NUMBER
	SKIPG	LPSA		;WAS THE ARGUMENT POSITIVE?
	MOVEI	A, 0		;NO, RETURN 0
	JRST	EXPXIT		;AND RETURN

EXP1:	MULI	LPSA,400	;SEPARATE FRACTION AND EXPONENT
	TSC	LPSA,LPSA	;GET A POSITIVE EXPONENT
	MUL	TEMP,E5		;FIXED POINT MULTIPLY BY LOG2(E)
	ASHC	TEMP,-242(LPSA)	;SEPARATE FRACTION AND INTEGER
	AOSG	TEMP		;ALGORITHM CALLS FOR MULT. BY 2
	AOS	TEMP		;ADJUST IF FRACTION WAS NEGATIVE
	HRRM	TEMP,B 		;SAVE FOR FUTURE SCALING
	JUMPG	USER,ASHH	;GO AHEAD IF ARG GREATER THAN 0
	TRNN	USER,377	;ALL THESE BITS 0?
	 JRST	 ASHH		;YES -- GO AHEAD
	ADDI	USER,200	;NO -- FIX UP
ASHH:	ASH	USER, -10	;MAKE ROOM FOR EXPONENT
	TLC	USER, 200000	;PUT 200 IN EXPONENT BITS
	FADB	USER, -1(P) 	;NORMALIZE, RESULTS TO USER AND E
	FMP	USER,USER	;FORM X↑2
	MOVE	A, E2		;GET FIRST CONSTANT
	FMP	A, USER		;E2*X↑2 IN A
	FAD	USER, E4	;ADD E4 TO RESULTS IN USER
	MOVE	LPSA, E3	;PICK UP E3
	FDV	LPSA,USER	;CALCULATE E3/(F↑2 + E4)
	FSB	A,LPSA		;E2*F↑2-E3(F↑2 + E4)**-1
	MOVE	TEMP,-1(P)  	;GET F AGAIN
	FSB	A, TEMP		;SUBTRACT FROM PARTIAL SUM
	FAD	A, E1		;ADD IN E1
	FDVM	TEMP, A		;DIVIDE BY F
	FAD	A, E6		;ADD 0.5
	FSC	A, (B)		;SCALE THE RESULTS
EXPXIT:	POP	P,B		;RESTORE AC
	SUB	P,X33		;ADJUST STACK
	JRST	@2(P)		;RETURN.

E1:	204476430062		;9.95459578
E2:	174433723400		;0.03465735903
E3:	212464770715		;617.97226953
E4:	207535527022		;87.417497202
E5:	270524354513		;LOG(E), BASE 2
E6:	0.5
E7:	207540071260		;88.028


;FLOATING POINT SINGLE PRECISION LOGARITHM FUNCTION
;LOG(ABSF(X)) IS CALCULATED BY THE SUBROUTINE, AND AN
;ARGUMENT OF ZERO IS RETURNED AS MINUS INFINITY. THE ALGORITHM IS

;LOGE(X) = (I + LOG2(F))*LOGE(2)
;WHERE X = (F/2)*2↑(I+1), AND LOG2(F) IS GIVEN BY
;LOG2(F) = C1*Z + C3*Z↑3 + C5*Z↑5 - 1/2
;AND Z = (F-SQRT(2))/(F+SQRT(2))

;THE CALLING SEQUENCE FOR THE ROUTINE IS AS FOLLOWS:
;	PUSH P,ARG
;	PUSHJ	P, LOG
;THE ANSWER IS RETURNED IN ACCUMULATOR A


HERE(LOG$)
.LOG:
	SKIPGE	-1(P)		;CHECK SIGN OF ARGUMENT.
	ERR <LOG: Negative argument -- real part returned>,1
	MOVM	LPSA,-1(P)   	;GET ABSF(A)
	JUMPE	LPSA, LZERO	;CHECK FOR ZERO ARGUMENT
	CAMN	LPSA, ONE	;CHECK FOR 1.0 ARGUMENT
	JRST	ZERANS		;IT IS 1.0 RETURN ZERO ANS.
	ASHC	LPSA, -33	;SEPARATE FRACTION FROM EXPONENT
	ADDI	LPSA, 211000	;FLOAT THE EXPONENT AND MULT. BY 2
	MOVSM	LPSA,USER	;NUMBER NOW IN CORRECT FL. FORMAT
	MOVSI	LPSA, 567377	;SET UP -401.0 IN LPSA
	FADM	LPSA,USER 	;SUBTRACT 401 FROM EXP.*2
	ASH	TEMP, -10	;SHIFT FRACTION FOR FLOATING
	TLC	TEMP, 200000	;FLOAT THE FRACTION PART
	FAD	TEMP, L1	;TEMP = TEMP-SQRT(2.0)/2.0
	MOVE	LPSA,TEMP	;PUT RESULTS IN LPSA
	FAD	LPSA, L2	;LPSA = LPSA+SQRT(2.0)
	FDV	TEMP,LPSA	;TEMP = TEMP/LPSA
	MOVEM	TEMP,A		;STORE NEW VARIABLE IN A
	FMP	TEMP,TEMP	;CALCULATE Z↑2
	MOVE	LPSA, L3	;PICK UP FIRST CONSTANT
	FMP	LPSA,TEMP	;MULTIPLY BY Z↑2
	FAD	LPSA, L4	;ADD IN NEXT CONSTANT
	FMP	LPSA,TEMP	;MULTIPLY BY Z↑2
	FAD	LPSA, L5	;ADD IN NEXT CONSTANT
	FMP	A,LPSA		;MULTIPLY BY Z
	FAD	A,USER		;ADD IN EXPONENT TO FORM LOG2(X)
	FMP	A, L7		;MULTIPLY TO FORM LOGE(X)
LOGXIT:	SUB	P,X22
	JRST	@2(P)

LZERO:	SKIPA	A, MIFI		;PICK UP MINUS INFINITY
ZERANS:	MOVEI	A,0		;MARG ANS ZERO
	JRST	LOGXIT		;AND RETURN

;CONSTANTS

ONE:	201400000000
L1:	577225754146		;-0.707106781187
L2:	201552023632		;1.414213562374
L3:	200462532521		;0.5989786496
L4:	200754213604		;0.9614706323
L5:	202561251002		;2.8853912903
L7:	200542710300		;0.69314718056
MIFI:	400000000001		;LARGEST NEGATIVE FLOATING NUMBER

ENDCOM (POW)


COMPIL(COD,<CODE,CALL>,<.SKIP.,CVSIX,X22,GOGTAB,X33>,<CODE, CALL>)